plot_ly(data = data_raw_0[(data_raw_0$`Sucide Rate` != 0) & (data_raw_0$`Urban population` != 0),],
x = ~ `Sucide Rate`,
y = ~ `Urban population`,
color = ~ continent,
marker = list(
size = 10
)
) %>%
layout(
title = 'Отношение Sucide Rate к Urban population',
yaxis = list(title = 'Urban population',
zeroline = FALSE), # Уберём выделения нулевых осей по y
xaxis = list(title = 'Sucide Rate',
zeroline = FALSE)) # Уберём выделения нулевых осей по y
data_raw_0 %>%
filter(continent == "Africa" | continent == "Americas") %>%
group_by(continent) %>%
get_summary_stats(`Life expectancy`, type = "mean_sd")%>%
ungroup()
## # A tibble: 2 × 5
## continent variable n mean sd
## <fct> <fct> <dbl> <dbl> <dbl>
## 1 Africa Life expectancy 52 66.1 5.97
## 2 Americas Life expectancy 38 78.6 3.68
data_raw_filtred <- data_raw_0 %>% filter(continent == "Africa" | continent == "Americas")
ggqqplot(data_raw_filtred[data_raw_filtred$`Life expectancy`,],
x = "Life expectancy", facet.by = "continent")
Так как данные распределены “не нормально” предлагаю использовать параметрический тест, например U-тест Манна-Уитни
library(rstatix)
# Создадим два датафрейма с значениями Life expectancy для разных континентов
group_af <- data_raw_0 %>% filter(continent == "Africa") %>% select(`Life expectancy`) %>% pull(`Life expectancy`)
group_am <- data_raw_0 %>% filter(continent == "Americas") %>% select(`Life expectancy`) %>% pull(`Life expectancy`)
# Проведем тест Манна-Уитни
p_valie_Utest <- wilcox.test(group_af, group_am)
p_valie_Utest
##
## Wilcoxon rank sum test with continuity correction
##
## data: group_af and group_am
## W = 107, p-value = 6.342e-13
## alternative hypothesis: true location shift is not equal to 0
wilcox_test_result <- wilcox.test(`Life expectancy` ~ continent, data = data_raw_filtred)
# Визуализация результатов
ggplot(data_raw_filtred, aes(x = continent, y = `Life expectancy`, fill = continent)) +
geom_boxplot() +
geom_dotplot(binaxis = "y", stackdir = "center", position = "dodge") +
stat_compare_means(test = "wilcox.test", label = "p.format", comparisons = list(c("Africa", "Americas"))) +
labs(title = "Comparison of Life Expectancy between Africa and Americas") +
theme_minimal()
# Выбираем все числовые колонки, кроме 'Year'
data_raw_0 <- data_raw_0 %>%
mutate(continent_num = case_when(continent == "Africa" ~ 1,
continent == "Americas" ~ 2,
TRUE ~ 0))
numerical_data <- data_raw_0 %>% select(where(is.numeric)) %>% filter(continent_num > 0)
numerical_data <- subset(numerical_data, select = -c(Year))
# Проводим корреляционный анализ
correlation_matrix <- cor(numerical_data)
# Визуализация корреляций - используем heatmap
melted_correlation <- melt(correlation_matrix)
ggplot(melted_correlation, aes(x = Var1, y = Var2, fill = value)) +
geom_tile(color = "white") +
scale_fill_gradient2(low = "blue", high = "red", mid = "white", midpoint = 0, limit = c(-1, 1), space = "Lab", name="Correlation") +
theme_minimal() +
theme(axis.text.x = element_text(angle = 45, vjust = 1, size = 8, hjust = 1)) +
labs(title = "Correlation Heatmap")
corrplot(correlation_matrix, method = 'ellipse', order = 'AOE', tl.cex = 0.5, tl.srt = 45)
numerical_data_scaled <- scale(numerical_data)
dist <- dist(numerical_data_scaled,
method = "euclidean"
)
as.matrix(dist)[1:6,1:6]
## 1 2 3 4 5 6
## 1 0.000000 6.074379 5.025993 5.170680 5.884635 3.816022
## 2 6.074379 0.000000 8.099954 7.772997 6.845781 6.879864
## 3 5.025993 8.099954 0.000000 5.590110 4.748777 4.197357
## 4 5.170680 7.772997 5.590110 0.000000 4.719991 3.555197
## 5 5.884635 6.845781 4.748777 4.719991 0.000000 4.741891
## 6 3.816022 6.879864 4.197357 3.555197 4.741891 0.000000
clear_hc <- hclust(d = dist, method = "ward.D2")
fviz_dend(clear_hc,
cex = 0.1)
## Warning: The `<scale>` argument of `guides()` cannot be `FALSE`. Use "none" instead as
## of ggplot2 3.3.4.
## ℹ The deprecated feature was likely used in the factoextra package.
## Please report the issue at <https://github.com/kassambara/factoextra/issues>.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
# Создаем heatmap с дендрограммой
pheatmap(numerical_data_scaled,
show_rownames = FALSE,
clustering_distance_rows = dist,
clustering_method = "ward.D2",
cutree_rows = 3,
cutree_cols = length(colnames(numerical_data_scaled)),
angle_col = 45,
main = "Dendrograms for clustering rows and columns with heatmap")
Делаем PCA:
full.pca <- prcomp(numerical_data_scaled,
scale = F)
summary(full.pca)
## Importance of components:
## PC1 PC2 PC3 PC4 PC5 PC6 PC7
## Standard deviation 2.7815 1.6834 1.40851 1.36879 1.10759 0.98271 0.84961
## Proportion of Variance 0.3868 0.1417 0.09919 0.09368 0.06134 0.04829 0.03609
## Cumulative Proportion 0.3868 0.5285 0.62774 0.72142 0.78275 0.83104 0.86713
## PC8 PC9 PC10 PC11 PC12 PC13 PC14
## Standard deviation 0.74926 0.70909 0.65772 0.55364 0.47697 0.43659 0.41107
## Proportion of Variance 0.02807 0.02514 0.02163 0.01533 0.01138 0.00953 0.00845
## Cumulative Proportion 0.89520 0.92034 0.94197 0.95730 0.96867 0.97820 0.98665
## PC15 PC16 PC17 PC18 PC19 PC20
## Standard deviation 0.3317 0.27505 0.22546 0.16374 0.06028 2.257e-16
## Proportion of Variance 0.0055 0.00378 0.00254 0.00134 0.00018 0.000e+00
## Cumulative Proportion 0.9921 0.99594 0.99848 0.99982 1.00000 1.000e+00
fviz_eig(full.pca, addlabels = T, ylim = c(0, 50), ncp = 19)
fviz_pca_var(full.pca,
select.var = list(contrib = 5), # Задаём число здесь
col.var = "contrib")
Наибольший вклад в “коррелированность” данных вносят Life expectancy, Infant Mortality, HepB3 Immunization, DPT Immunization, Measles Immunization
fviz_contrib(full.pca, choice = "var", axes = 1, top = 30)
fviz_contrib(full.pca, choice = "var", axes = 2, top = 24)
fviz_contrib(full.pca, choice = "var", axes = 3, top = 24)
fviz_contrib(full.pca, choice = "var", axes = 4, top = 24)
fviz_contrib(full.pca, choice = "var", axes = 5, top = 24)
Рассмотрим Dim-1. В данных нет ярковыраженных переменных, которые бы брали на себя наибольших “вклад” для первой главной компоненты однако можно выделить группу переменных, которые отвечают за приблизительно 40% такого вклада:
Life expectancy
Basic sanitation services
Infant Mortality
Clean fuels and cooking technologies
# Визуализируем с группировкой по continent (для этого переменную нужно сделать фактором)
ggbiplot(full.pca,
scale=0,
groups = as.factor(numerical_data$continent_num),
ellipse = T,
alpha = 0.2) +
theme_minimal()
Значащие группы 1 и 2 – это Африка и Америки соответственно.
Если мы доверяем построиным линиям то можно сформулировать следубщие гипотрезы:
Страны с большим Sucide Rate более вероятно будут относиться к континенту африка
Страны с больней Life expectancy более вероятно относятся к континантам Америки
Rural population больше в странах Африки, Urban population больше в странах относящихся к Америке
Показатель Clean fuels and cooking technologies и Basic sanitation services более выражены в странах Америки
и другие
#install.packages(c("irlba", "Matrix"))
library(tidymodels)
library(embed)
umap_prep <- recipe(~., data = numerical_data) %>% # "техническая" строка, нужная для работы фреймворка tidymodels
step_normalize(all_predictors()) %>% # нормируем все колонки
step_umap(all_predictors()) %>% # проводим в UMAP. Используем стандартные настройки. Чтобы менять ключевой параметр (neighbors), нужно больше погружаться в машинное обучение
prep() %>% # "техническая" строка, нужная для работы фреймворка tidymodels. Мы выполняем все степы выше
juice() # Финальная строка - приводим результаты UMAP к стандартизированному датасету
umap_prep %>%
ggplot(aes(UMAP1, UMAP2)) + # # можно добавить раскраску
geom_point(aes(color = as.factor(numerical_data$continent_num),
shape = numerical_data$diabetes_ch),
alpha = 0.7, size = 2) +
labs(color = NULL)+
ggtitle("Отображение данных в зависимости от Континента")
Сложно интерпритировать UMAP график. Однако если сравнивать его с PCA графиком, то можно заметить что:
PCA: Основан на линейных преобразованиях и стремится сохранить максимальную дисперсию данных в новых компонентах (главных компонентах).
UMAP: Основан на нелинейных преобразованиях, стремится сохранить локальные структуры данных, предоставляя более низкоразмерное представление, сохраняя при этом глобальные структуры.
PCA: Сохраняет глобальные расстояния, что может привести к сжатию или искажению локальных структур данных.
UMAP: Стремится сохранять локальные расстояния, что делает его лучшим для визуализации кластеров и плотных областей в данных.
numerical_data_scaled_drop_1_5 <- numerical_data_scaled[, -c(1:5)]
numerical_data_scaled_drop_6_10 <- numerical_data_scaled[, -c(6:10)]
numerical_data_scaled_drop_11_15 <- numerical_data_scaled[, -c(11:15)]
full_drop_5.pca <- prcomp(numerical_data_scaled_drop_1_5,
scale = F)
summary(full_drop_5.pca)
## Importance of components:
## PC1 PC2 PC3 PC4 PC5 PC6 PC7
## Standard deviation 2.4452 1.6408 1.2933 1.04801 0.86905 0.84416 0.77616
## Proportion of Variance 0.3986 0.1795 0.1115 0.07322 0.05035 0.04751 0.04016
## Cumulative Proportion 0.3986 0.5781 0.6896 0.76281 0.81316 0.86067 0.90083
## PC8 PC9 PC10 PC11 PC12 PC13 PC14
## Standard deviation 0.68995 0.58551 0.50801 0.45351 0.31760 0.27616 0.16689
## Proportion of Variance 0.03174 0.02285 0.01721 0.01371 0.00672 0.00508 0.00186
## Cumulative Proportion 0.93256 0.95542 0.97262 0.98633 0.99306 0.99814 1.00000
## PC15
## Standard deviation 2.294e-16
## Proportion of Variance 0.000e+00
## Cumulative Proportion 1.000e+00
full_drop_10.pca <- prcomp(numerical_data_scaled_drop_6_10,
scale = F)
summary(full_drop_10.pca)
## Importance of components:
## PC1 PC2 PC3 PC4 PC5 PC6 PC7
## Standard deviation 2.3863 1.4573 1.3644 1.2733 1.01672 0.91987 0.80672
## Proportion of Variance 0.3796 0.1416 0.1241 0.1081 0.06891 0.05641 0.04339
## Cumulative Proportion 0.3796 0.5212 0.6453 0.7534 0.82231 0.87872 0.92211
## PC8 PC9 PC10 PC11 PC12 PC13 PC14
## Standard deviation 0.65115 0.48608 0.43437 0.40561 0.29331 0.25535 0.06080
## Proportion of Variance 0.02827 0.01575 0.01258 0.01097 0.00574 0.00435 0.00025
## Cumulative Proportion 0.95037 0.96612 0.97870 0.98967 0.99541 0.99975 1.00000
## PC15
## Standard deviation 2.879e-16
## Proportion of Variance 0.000e+00
## Cumulative Proportion 1.000e+00
full_drop_15.pca <- prcomp(numerical_data_scaled_drop_11_15,
scale = F)
summary(full_drop_15.pca)
## Importance of components:
## PC1 PC2 PC3 PC4 PC5 PC6 PC7
## Standard deviation 2.456 1.464 1.3556 1.15139 1.02849 0.82952 0.72496
## Proportion of Variance 0.402 0.143 0.1225 0.08838 0.07052 0.04587 0.03504
## Cumulative Proportion 0.402 0.545 0.6675 0.75584 0.82636 0.87223 0.90727
## PC8 PC9 PC10 PC11 PC12 PC13 PC14
## Standard deviation 0.65098 0.60165 0.45755 0.40791 0.39998 0.25634 0.06107
## Proportion of Variance 0.02825 0.02413 0.01396 0.01109 0.01067 0.00438 0.00025
## Cumulative Proportion 0.93552 0.95966 0.97361 0.98471 0.99537 0.99975 1.00000
## PC15
## Standard deviation 2.384e-16
## Proportion of Variance 0.000e+00
## Cumulative Proportion 1.000e+00
Действительно Cumulative Proportion отличается в зависимости от выбранных переменных. Отличия не большие +- 5 процентных пунктов в моем случае.
Давайте самостоятельно увидим, что снижение размерности – это группа методов, славящаяся своей неустойчивостью. Создайте две дамми-колонки о том: (1) принадлежит ли страна к африканскому континенту, (2) Океании. Проведите PCA вместе с ними, постройте биплоты. Проинтерпрейтируйте результат. Объясните, почему добавление дамми-колонок не совсем корректно в случае PCA нашего типа.
data_raw_dummy <- data_raw_0 %>%
mutate(is_africa = ifelse(continent== "Africa", 1, 0),
is_americas = ifelse(continent == "Americas", 1, 0)) %>%
filter(continent_num > 0)
data_raw_dummy <- subset(data_raw_dummy, select = -c(Year)) %>%
select(where(is.numeric))
data_raw_dummy_scaled <- scale(data_raw_dummy)
dist <- dist(data_raw_dummy_scaled,
method = "euclidean"
)
as.matrix(dist)[1:6,1:6]
## 1 2 3 4 5 6
## 1 0.000000 6.074379 5.776495 5.902815 6.537294 4.761230
## 2 6.074379 0.000000 8.585834 8.278090 7.414311 7.445792
## 3 5.776495 8.585834 0.000000 5.590110 4.748777 4.197357
## 4 5.902815 8.278090 5.590110 0.000000 4.719991 3.555197
## 5 6.537294 7.414311 4.748777 4.719991 0.000000 4.741891
## 6 4.761230 7.445792 4.197357 3.555197 4.741891 0.000000
dummy.pca <- prcomp(data_raw_dummy_scaled,
scale = F)
summary(dummy.pca)
## Importance of components:
## PC1 PC2 PC3 PC4 PC5 PC6 PC7
## Standard deviation 3.0154 1.6938 1.41620 1.38738 1.22256 0.98898 0.88165
## Proportion of Variance 0.4133 0.1304 0.09116 0.08749 0.06794 0.04446 0.03533
## Cumulative Proportion 0.4133 0.5437 0.63488 0.72237 0.79031 0.83477 0.87010
## PC8 PC9 PC10 PC11 PC12 PC13 PC14
## Standard deviation 0.79926 0.71470 0.65792 0.57864 0.55150 0.43942 0.41479
## Proportion of Variance 0.02904 0.02322 0.01968 0.01522 0.01383 0.00878 0.00782
## Cumulative Proportion 0.89914 0.92236 0.94203 0.95725 0.97108 0.97985 0.98767
## PC15 PC16 PC17 PC18 PC19 PC20
## Standard deviation 0.3349 0.27870 0.22559 0.16379 0.06028 2.534e-16
## Proportion of Variance 0.0051 0.00353 0.00231 0.00122 0.00017 0.000e+00
## Cumulative Proportion 0.9928 0.99630 0.99862 0.99983 1.00000 1.000e+00
## PC21 PC22
## Standard deviation 1.47e-16 3.329e-32
## Proportion of Variance 0.00e+00 0.000e+00
## Cumulative Proportion 1.00e+00 1.000e+00
fviz_eig(dummy.pca, addlabels = T, ylim = c(0, 50), ncp = 19)
первая и последующие переменные стали “объяснять” меньший процент корреляции данных
fviz_pca_var(dummy.pca,
select.var = list(contrib = 5), # Задаём число здесь
col.var = "contrib")
первые 5 переменных не изменились но направление изменилось
# Визуализируем с группировкой по continent (для этого переменную нужно сделать фактором)
ggbiplot(dummy.pca,
scale=0,
groups = as.factor(data_raw_dummy$continent_num),
ellipse = T,
alpha = 0.2) +
theme_minimal()
Добавив бинарные дамми колонки которые характеризуют наблюдения в зависимости от пренадлежности к материку: Африка или Америка(Я не стал брать Океанию для консистенси с прошлыми графиками) похоже что таким образом мы “отделили” наблюдения одного материка от другого, но таким образом потеряли визуальное подтверждение “схожести” некоторых наблюдений (пересечение овалов на предыдущем ggbiplot)